The data has been collected during the study of . The dataset is available in open-access at …
library(tidyverse)
library(ggside)
library(easystats)
library(patchwork)
library(mirt)
df <- read.csv("https://raw.githubusercontent.com/RealityBending/IllusionGameReliability/main/data/preprocessed_questionnaires.csv") |>
# df <- read.csv("../IllusionGameReliability/data/preprocessed_questionnaires.csv") |>
select(Sex, Age, starts_with("Item_PHQ4")) |>
filter(!is.na(Item_PHQ4_Anxiety_1))
names(df) <- str_remove_all(names(df), "Item_PHQ4_")
report::report_participants(df)
## [1] "485 participants (Mean age = 30.1, SD = 10.1, range: [18, 73]; Sex: 50.3% females, 49.7% males, 0.0% other)"
add_labels <- function(x) {
x <- case_when(
x == 0 ~ "Not at all",
x == 1 ~ "Once or twice",
x == 2 ~ "Several days",
x == 3 ~ "More than half the days",
TRUE ~ "Nearly every day"
)
fct_relevel(x, c("Not at all", "Once or twice", "Several days", "More than half the days", "Nearly every day"))
}
df <- select(df, -Sex, -Age)
data <- df |>
pivot_longer(everything(), names_to = "Item", values_to = "Answer") |>
group_by(Item, Answer) |>
summarise(n = n() / nrow(df)) |>
mutate(
Facet = str_split_fixed(Item, "_", 2)[, 1],
# Item = str_split_fixed(Item, "_", 2)[,2],
Answer = add_labels(Answer),
Item = case_when(
Item == "Anxiety_1" ~ "A1: Feeling nervous, anxious or on edge",
Item == "Anxiety_2" ~ "A2: Not being able to stop or control worrying",
Item == "Depression_3" ~ "D1: Little interest or pleasure in doing things",
TRUE ~ "D2: Feeling down, depressed, or hopeless"
))
p1 <- data |>
ggplot(aes(x = Answer, fill=Facet)) +
geom_bar(aes(y = n), stat = "identity") +
scale_y_continuous(labels=scales::percent, expand = c(0, 0)) +
scale_fill_manual(values = c("Depression" = "#3F51B5", "Anxiety" = "#E91E63"), guide = "none") +
labs(y = "Proportion of Answers", title = "A. Proportion of answers per item") +
theme_modern(axis.title.space = 10) +
theme(axis.text.x = element_text(angle = 45, hjust=1),
axis.title.x = element_blank(),
plot.title = element_text(face = "bold"),
strip.background=element_rect(fill="#EEEEEE", colour="white")) +
facet_wrap(~Item)
p1
The “Once or twice” answer was selected in 29.85%
p1b <- rbind(
data.frame(Scores = paste0(df$Anxiety_1, "_", df$Anxiety_2),
Facet = "Anxiety"),
data.frame(Scores = paste0(df$Depression_3, "_", df$Depression_4),
Facet = "Depression")
) |>
group_by(Facet, Scores) |>
summarize(n = n() / nrow(df)) |>
separate(Scores, into = c("Q1", "Q2")) |>
mutate(Label = ifelse(Q1 == Q2, format_percent(n), "")) |>
mutate(Q1 = add_labels(as.numeric(Q1)),
Q2 = add_labels(as.numeric(Q2))) |>
ggplot(aes(x = Q1, y = Q2)) +
geom_tile(aes(fill = n)) +
geom_text(aes(label = Label)) +
scale_y_discrete(expand = c(0, 0)) +
scale_x_discrete(expand = c(0, 0)) +
scale_fill_gradientn(colors=c("white", "#2196F3", "#3F51B5", "#673AB7"), labels = scales::percent) +
labs(fill = "% Participants", x = "Item 1", y = "Item 2", title = "Joint prevalence of responses") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust=1),
panel.grid.major = element_blank(),
plot.title = element_text(face = "bold"),
strip.background=element_rect(fill="#EEEEEE", colour="white"),
strip.text = element_text(face = "bold")) +
facet_wrap(~Facet)
p1b
df_anx <- select(df, contains("Anxiety"))
performance::cronbachs_alpha(df_anx)
## [1] 0.903
m_anxiety <- mirt(df_anx, model = 1, itemtype = 'graded', SE = TRUE)
##
Iteration: 1, Log-Lik: -1434.510, Max-Change: 0.49722
Iteration: 2, Log-Lik: -1377.689, Max-Change: 0.62182
Iteration: 3, Log-Lik: -1322.367, Max-Change: 0.72396
Iteration: 4, Log-Lik: -1283.817, Max-Change: 0.66639
Iteration: 5, Log-Lik: -1263.174, Max-Change: 0.50226
Iteration: 6, Log-Lik: -1253.813, Max-Change: 0.34052
Iteration: 7, Log-Lik: -1249.772, Max-Change: 0.23422
Iteration: 8, Log-Lik: -1247.939, Max-Change: 0.16157
Iteration: 9, Log-Lik: -1247.027, Max-Change: 0.11472
Iteration: 10, Log-Lik: -1246.280, Max-Change: 0.06322
Iteration: 11, Log-Lik: -1246.106, Max-Change: 0.05105
Iteration: 12, Log-Lik: -1246.000, Max-Change: 0.04207
Iteration: 13, Log-Lik: -1245.837, Max-Change: 0.01752
Iteration: 14, Log-Lik: -1245.829, Max-Change: 0.01703
Iteration: 15, Log-Lik: -1245.823, Max-Change: 0.01663
Iteration: 16, Log-Lik: -1245.792, Max-Change: 0.01528
Iteration: 17, Log-Lik: -1245.788, Max-Change: 0.01454
Iteration: 18, Log-Lik: -1245.783, Max-Change: 0.01500
Iteration: 19, Log-Lik: -1245.760, Max-Change: 0.01522
Iteration: 20, Log-Lik: -1245.755, Max-Change: 0.01616
Iteration: 21, Log-Lik: -1245.750, Max-Change: 0.01580
Iteration: 22, Log-Lik: -1245.726, Max-Change: 0.01663
Iteration: 23, Log-Lik: -1245.722, Max-Change: 0.01803
Iteration: 24, Log-Lik: -1245.718, Max-Change: 0.01645
Iteration: 25, Log-Lik: -1245.694, Max-Change: 0.02301
Iteration: 26, Log-Lik: -1245.689, Max-Change: 0.01862
Iteration: 27, Log-Lik: -1245.684, Max-Change: 0.01712
Iteration: 28, Log-Lik: -1245.660, Max-Change: 0.02176
Iteration: 29, Log-Lik: -1245.656, Max-Change: 0.02120
Iteration: 30, Log-Lik: -1245.651, Max-Change: 0.02021
Iteration: 31, Log-Lik: -1245.626, Max-Change: 0.01921
Iteration: 32, Log-Lik: -1245.622, Max-Change: 0.02052
Iteration: 33, Log-Lik: -1245.618, Max-Change: 0.01939
Iteration: 34, Log-Lik: -1245.594, Max-Change: 0.02217
Iteration: 35, Log-Lik: -1245.589, Max-Change: 0.02285
Iteration: 36, Log-Lik: -1245.585, Max-Change: 0.02256
Iteration: 37, Log-Lik: -1245.561, Max-Change: 0.01765
Iteration: 38, Log-Lik: -1245.557, Max-Change: 0.01948
Iteration: 39, Log-Lik: -1245.553, Max-Change: 0.02142
Iteration: 40, Log-Lik: -1245.529, Max-Change: 0.02326
Iteration: 41, Log-Lik: -1245.525, Max-Change: 0.02321
Iteration: 42, Log-Lik: -1245.521, Max-Change: 0.02237
Iteration: 43, Log-Lik: -1245.499, Max-Change: 0.02324
Iteration: 44, Log-Lik: -1245.495, Max-Change: 0.02257
Iteration: 45, Log-Lik: -1245.491, Max-Change: 0.02463
Iteration: 46, Log-Lik: -1245.468, Max-Change: 0.02322
Iteration: 47, Log-Lik: -1245.464, Max-Change: 0.02381
Iteration: 48, Log-Lik: -1245.461, Max-Change: 0.02434
Iteration: 49, Log-Lik: -1245.440, Max-Change: 0.02352
Iteration: 50, Log-Lik: -1245.436, Max-Change: 0.02455
Iteration: 51, Log-Lik: -1245.432, Max-Change: 0.02511
Iteration: 52, Log-Lik: -1245.411, Max-Change: 0.02529
Iteration: 53, Log-Lik: -1245.408, Max-Change: 0.02560
Iteration: 54, Log-Lik: -1245.405, Max-Change: 0.02584
Iteration: 55, Log-Lik: -1245.385, Max-Change: 0.02519
Iteration: 56, Log-Lik: -1245.381, Max-Change: 0.02795
Iteration: 57, Log-Lik: -1245.378, Max-Change: 0.02268
Iteration: 58, Log-Lik: -1245.365, Max-Change: 0.02836
Iteration: 59, Log-Lik: -1245.361, Max-Change: 0.02557
Iteration: 60, Log-Lik: -1245.358, Max-Change: 0.02714
Iteration: 61, Log-Lik: -1245.340, Max-Change: 0.02686
Iteration: 62, Log-Lik: -1245.337, Max-Change: 0.02687
Iteration: 63, Log-Lik: -1245.334, Max-Change: 0.02683
Iteration: 64, Log-Lik: -1245.316, Max-Change: 0.02909
Iteration: 65, Log-Lik: -1245.313, Max-Change: 0.02611
Iteration: 66, Log-Lik: -1245.311, Max-Change: 0.02433
Iteration: 67, Log-Lik: -1245.296, Max-Change: 0.02794
Iteration: 68, Log-Lik: -1245.293, Max-Change: 0.02868
Iteration: 69, Log-Lik: -1245.290, Max-Change: 0.02824
Iteration: 70, Log-Lik: -1245.274, Max-Change: 0.03027
Iteration: 71, Log-Lik: -1245.272, Max-Change: 0.03022
Iteration: 72, Log-Lik: -1245.269, Max-Change: 0.02880
Iteration: 73, Log-Lik: -1245.253, Max-Change: 0.01696
Iteration: 74, Log-Lik: -1245.252, Max-Change: 0.03171
Iteration: 75, Log-Lik: -1245.249, Max-Change: 0.02984
Iteration: 76, Log-Lik: -1245.238, Max-Change: 0.02126
Iteration: 77, Log-Lik: -1245.234, Max-Change: 0.02745
Iteration: 78, Log-Lik: -1245.231, Max-Change: 0.00321
Iteration: 79, Log-Lik: -1245.231, Max-Change: 0.02698
Iteration: 80, Log-Lik: -1245.229, Max-Change: 0.02876
Iteration: 81, Log-Lik: -1245.226, Max-Change: 0.02913
Iteration: 82, Log-Lik: -1245.213, Max-Change: 0.02962
Iteration: 83, Log-Lik: -1245.211, Max-Change: 0.02906
Iteration: 84, Log-Lik: -1245.209, Max-Change: 0.02921
Iteration: 85, Log-Lik: -1245.196, Max-Change: 0.03105
Iteration: 86, Log-Lik: -1245.194, Max-Change: 0.02980
Iteration: 87, Log-Lik: -1245.192, Max-Change: 0.03064
Iteration: 88, Log-Lik: -1245.179, Max-Change: 0.02993
Iteration: 89, Log-Lik: -1245.177, Max-Change: 0.03042
Iteration: 90, Log-Lik: -1245.175, Max-Change: 0.03243
Iteration: 91, Log-Lik: -1245.164, Max-Change: 0.02865
Iteration: 92, Log-Lik: -1245.161, Max-Change: 0.03030
Iteration: 93, Log-Lik: -1245.159, Max-Change: 0.03026
Iteration: 94, Log-Lik: -1245.149, Max-Change: 0.02627
Iteration: 95, Log-Lik: -1245.147, Max-Change: 0.00253
Iteration: 96, Log-Lik: -1245.147, Max-Change: 0.03315
Iteration: 97, Log-Lik: -1245.145, Max-Change: 0.02951
Iteration: 98, Log-Lik: -1245.144, Max-Change: 0.03432
Iteration: 99, Log-Lik: -1245.142, Max-Change: 0.03095
Iteration: 100, Log-Lik: -1245.134, Max-Change: 0.02948
Iteration: 101, Log-Lik: -1245.132, Max-Change: 0.03010
Iteration: 102, Log-Lik: -1245.131, Max-Change: 0.03003
Iteration: 103, Log-Lik: -1245.122, Max-Change: 0.02659
Iteration: 104, Log-Lik: -1245.120, Max-Change: 0.02996
Iteration: 105, Log-Lik: -1245.118, Max-Change: 0.03053
Iteration: 106, Log-Lik: -1245.110, Max-Change: 0.03140
Iteration: 107, Log-Lik: -1245.108, Max-Change: 0.03240
Iteration: 108, Log-Lik: -1245.107, Max-Change: 0.03022
Iteration: 109, Log-Lik: -1245.099, Max-Change: 0.03116
Iteration: 110, Log-Lik: -1245.097, Max-Change: 0.02909
Iteration: 111, Log-Lik: -1245.096, Max-Change: 0.03552
Iteration: 112, Log-Lik: -1245.091, Max-Change: 0.02542
Iteration: 113, Log-Lik: -1245.087, Max-Change: 0.02774
Iteration: 114, Log-Lik: -1245.085, Max-Change: 0.00788
Iteration: 115, Log-Lik: -1245.085, Max-Change: 0.03110
Iteration: 116, Log-Lik: -1245.084, Max-Change: 0.03182
Iteration: 117, Log-Lik: -1245.082, Max-Change: 0.03100
Iteration: 118, Log-Lik: -1245.075, Max-Change: 0.02973
Iteration: 119, Log-Lik: -1245.074, Max-Change: 0.00171
Iteration: 120, Log-Lik: -1245.074, Max-Change: 0.03305
Iteration: 121, Log-Lik: -1245.073, Max-Change: 0.03095
Iteration: 122, Log-Lik: -1245.072, Max-Change: 0.03123
Iteration: 123, Log-Lik: -1245.070, Max-Change: 0.02556
Iteration: 124, Log-Lik: -1245.066, Max-Change: 0.03212
Iteration: 125, Log-Lik: -1245.064, Max-Change: 0.03086
Iteration: 126, Log-Lik: -1245.063, Max-Change: 0.02970
Iteration: 127, Log-Lik: -1245.057, Max-Change: 0.02988
Iteration: 128, Log-Lik: -1245.056, Max-Change: 0.02452
Iteration: 129, Log-Lik: -1245.055, Max-Change: 0.02659
Iteration: 130, Log-Lik: -1245.050, Max-Change: 0.03050
Iteration: 131, Log-Lik: -1245.049, Max-Change: 0.00127
Iteration: 132, Log-Lik: -1245.049, Max-Change: 0.02366
Iteration: 133, Log-Lik: -1245.048, Max-Change: 0.03117
Iteration: 134, Log-Lik: -1245.047, Max-Change: 0.03019
Iteration: 135, Log-Lik: -1245.046, Max-Change: 0.02695
Iteration: 136, Log-Lik: -1245.042, Max-Change: 0.02767
Iteration: 137, Log-Lik: -1245.041, Max-Change: 0.02447
Iteration: 138, Log-Lik: -1245.040, Max-Change: 0.02816
Iteration: 139, Log-Lik: -1245.036, Max-Change: 0.02278
Iteration: 140, Log-Lik: -1245.035, Max-Change: 0.00219
Iteration: 141, Log-Lik: -1245.035, Max-Change: 0.03010
Iteration: 142, Log-Lik: -1245.034, Max-Change: 0.02924
Iteration: 143, Log-Lik: -1245.033, Max-Change: 0.02585
Iteration: 144, Log-Lik: -1245.032, Max-Change: 0.02333
Iteration: 145, Log-Lik: -1245.029, Max-Change: 0.03156
Iteration: 146, Log-Lik: -1245.028, Max-Change: 0.02978
Iteration: 147, Log-Lik: -1245.028, Max-Change: 0.00251
Iteration: 148, Log-Lik: -1245.028, Max-Change: 0.00250
Iteration: 149, Log-Lik: -1245.028, Max-Change: 0.03614
Iteration: 150, Log-Lik: -1245.027, Max-Change: 0.03004
Iteration: 151, Log-Lik: -1245.023, Max-Change: 0.00262
Iteration: 152, Log-Lik: -1245.023, Max-Change: 0.02886
Iteration: 153, Log-Lik: -1245.022, Max-Change: 0.02489
Iteration: 154, Log-Lik: -1245.019, Max-Change: 0.03097
Iteration: 155, Log-Lik: -1245.018, Max-Change: 0.00207
Iteration: 156, Log-Lik: -1245.018, Max-Change: 0.00074
Iteration: 157, Log-Lik: -1245.018, Max-Change: 0.00062
Iteration: 158, Log-Lik: -1245.018, Max-Change: 0.00066
Iteration: 159, Log-Lik: -1245.018, Max-Change: 0.00015
Iteration: 160, Log-Lik: -1245.018, Max-Change: 0.00013
Iteration: 161, Log-Lik: -1245.018, Max-Change: 0.00049
Iteration: 162, Log-Lik: -1245.018, Max-Change: 0.00059
Iteration: 163, Log-Lik: -1245.018, Max-Change: 0.00019
Iteration: 164, Log-Lik: -1245.018, Max-Change: 0.00040
Iteration: 165, Log-Lik: -1245.018, Max-Change: 0.00014
Iteration: 166, Log-Lik: -1245.018, Max-Change: 0.00012
Iteration: 167, Log-Lik: -1245.018, Max-Change: 0.00048
Iteration: 168, Log-Lik: -1245.018, Max-Change: 0.00010
##
## Calculating information matrix...
m_anxiety
##
## Call:
## mirt(data = df_anx, model = 1, itemtype = "graded", SE = TRUE)
##
## Full-information item factor analysis with 1 factor(s).
## Converged within 1e-04 tolerance after 168 EM iterations.
## mirt version: 1.37.1
## M-step optimizer: BFGS
## EM acceleration: Ramsay
## Number of rectangular quadrature: 61
## Latent density type: Gaussian
##
## Information matrix estimated with method: Oakes
## Second-order test: model is a possible local maximum
## Condition number of information matrix = 30371
##
## Log-likelihood = -1245
## Estimated parameters: 10
## AIC = 2510
## BIC = 2552; SABIC = 2520
## G2 (14) = 20.6, p = 0.113
## RMSEA = 0.031, CFI = NaN, TLI = NaN
# Factor loadings
summary(m_anxiety)
## F1 h2
## Anxiety_1 0.895 0.802
## Anxiety_2 0.991 0.982
##
## SS loadings: 1.78
## Proportion Var: 0.892
##
## Factor correlations:
##
## F1
## F1 1
# Alpha
a <- coef(m_anxiety)
a <- data.frame(
Item = c("Anxiety 1", "Anxiety 2"),
Discrimination = c(paste0(format_value(a$Anxiety_1[1, 1]), ", ", format_ci(a$Anxiety_1[2, 1], a$Anxiety_1[3, 1])), paste0(format_value(a$Anxiety_2[1, 1]), ", ", format_ci(a$Anxiety_2[2, 1], a$Anxiety_2[3, 1]))
))
knitr::kable(a)
| Item | Discrimination |
|---|---|
| Anxiety 1 | 3.42, 95% CI [2.66, 4.18] |
| Anxiety 2 | 12.55, 95% CI [-9.55, 34.65] |
# Plots
plot(m_anxiety, type = 'trace', theta_lim = c(-3, 3))
plot(m_anxiety, type = 'infotrace', theta_lim = c(-3, 3))
plot(m_anxiety, type = 'score', theta_lim = c(-3, 3))
plot(m_anxiety, type = 'infoSE', theta_lim = c(-3, 3))
plot(m_anxiety, type = 'rxx', theta_lim = c(-3, 3))
Typically, an optimally informative item will have a large location and broad category coverage (as indicated by location parameters) over theta.
crc <- function(mod, data) {
Theta <- matrix(seq(-3.5, 3.5, by = .01))
rez <- data.frame()
for(i in 1:2) {
dat <- as.data.frame(probtrace(extract.item(mod, i), Theta))
dat$Theta <- Theta[,1]
dat$Information <- normalize(iteminfo(extract.item(mod, i), Theta))
dat <- pivot_longer(dat, -one_of(c("Theta", "Information")), names_to = "Answer", values_to = "Probability")
dat$Item <- names(data)[i]
rez <- rbind(rez, dat)
}
rez <- rez |>
mutate(
Answer = case_when(
Answer == "P.1" ~ "Not at all",
Answer == "P.2" ~ "Once or twice",
Answer == "P.3" ~ "Several days",
Answer == "P.4" ~ "More than half the days",
TRUE ~ "Nearly every day"
),
Answer = fct_relevel(Answer, c("Not at all", "Once or twice", "Several days", "More than half the days", "Nearly every day")),
Item = case_when(
Item == "Anxiety_1" ~ "A1: Feeling nervous, anxious or on edge",
Item == "Anxiety_2" ~ "A2: Not being able to stop or control worrying",
Item == "Depression_3" ~ "D1: Little interest or pleasure in doing things",
TRUE ~ "D2: Feeling down, depressed, or hopeless"
))
# Get close to 0.95
i1 <- rez[rez$Item == unique(rez$Item)[1],]
i2 <- rez[rez$Item == unique(rez$Item)[2],]
minmax <- rbind(
rbind(
i1[i1$Answer == "Not at all", ][which.min(abs(i1[i1$Answer == "Not at all", ]$Probability - 0.95)), ],
i1[i1$Answer == "Nearly every day", ][which.min(abs(i1[i1$Answer == "Nearly every day", ]$Probability - 0.95)), ]
),
rbind(
i2[i2$Answer == "Not at all", ][which.min(abs(i2[i2$Answer == "Not at all", ]$Probability - 0.95)), ],
i2[i2$Answer == "Nearly every day", ][which.min(abs(i2[i2$Answer == "Nearly every day", ]$Probability - 0.95)), ]
)
)
p <- rez |>
ggplot(aes(x = Theta, y = Probability)) +
geom_line(aes(y = Information), linetype = "dotted", color = "grey") +
geom_line(aes(color = Answer), size=1) +
# geom_vline(data = minmax, aes(xintercept = Theta), linetype = "dashed") +
scale_y_continuous(labels=scales::percent, expand = c(0, 0.005)) +
scale_color_flat_d("rainbow") +
facet_grid(~Item) +
theme_modern(axis.title.space = 5) +
theme(strip.background=element_rect(fill="#EEEEEE", colour="white"))
list(p = p, rez = rez, minmax=minmax)
}
out <- crc(m_anxiety, df_anx)
p2a <- out$p + labs(x = expression(Anxiety~(theta)))
p2a
normalize_scores <- function(data, minmax) {
minmax <- minmax[minmax$Theta %in% range(minmax$Theta), ]
item <- data.frame()
scores <- data.frame()
for(i in unique(data$Item)) {
item <- data[data$Item == i, ] |>
filter(Theta >= min(minmax$Theta) & Theta <= max(minmax$Theta)) |>
mutate(Theta_n = normalize(Theta)) |>
rbind(item)
scores <- item[item$Item == i, ] |>
group_by(Answer) |>
filter(Probability == max(Probability)) |>
ungroup() |>
mutate(label = insight::format_value(Theta_n)) |>
rbind(scores)
}
item |>
ggplot(aes(x = Theta, y = Probability, color = Answer)) +
geom_line(size = 1) +
geom_segment(data=scores, aes(x=Theta, xend = Theta, y=0, yend=Probability, color = Answer), linetype = "dashed") +
geom_label(data=scores, aes(x=Theta, y=Probability, label=label), show.legend = FALSE) +
scale_y_continuous(labels=scales::percent, expand = c(0, 0.01)) +
scale_color_viridis_d(option = "D") +
facet_grid(~Item) +
theme_modern(axis.title.space = 5) +
theme(strip.background=element_rect(fill="#EEEEEE", colour="white"))
}
p3a <- normalize_scores(data=out$rez, out$minmax) + labs(x = expression(Anxiety~(theta)))
p3a
df_dep <- select(df, contains("Depression"))
performance::cronbachs_alpha(df_dep)
## [1] 0.841
m_depression <- mirt(df_dep, model = 1, itemtype = 'graded', SE = TRUE)
##
Iteration: 1, Log-Lik: -1454.429, Max-Change: 0.43021
Iteration: 2, Log-Lik: -1412.776, Max-Change: 0.48366
Iteration: 3, Log-Lik: -1375.753, Max-Change: 0.50038
Iteration: 4, Log-Lik: -1351.861, Max-Change: 0.41059
Iteration: 5, Log-Lik: -1339.730, Max-Change: 0.29809
Iteration: 6, Log-Lik: -1334.341, Max-Change: 0.20880
Iteration: 7, Log-Lik: -1332.019, Max-Change: 0.14321
Iteration: 8, Log-Lik: -1330.986, Max-Change: 0.09907
Iteration: 9, Log-Lik: -1330.503, Max-Change: 0.07020
Iteration: 10, Log-Lik: -1330.136, Max-Change: 0.03563
Iteration: 11, Log-Lik: -1330.079, Max-Change: 0.02845
Iteration: 12, Log-Lik: -1330.048, Max-Change: 0.02329
Iteration: 13, Log-Lik: -1330.002, Max-Change: 0.01230
Iteration: 14, Log-Lik: -1329.997, Max-Change: 0.01214
Iteration: 15, Log-Lik: -1329.991, Max-Change: 0.01163
Iteration: 16, Log-Lik: -1329.957, Max-Change: 0.01344
Iteration: 17, Log-Lik: -1329.950, Max-Change: 0.01509
Iteration: 18, Log-Lik: -1329.943, Max-Change: 0.01385
Iteration: 19, Log-Lik: -1329.900, Max-Change: 0.02264
Iteration: 20, Log-Lik: -1329.890, Max-Change: 0.01837
Iteration: 21, Log-Lik: -1329.882, Max-Change: 0.01814
Iteration: 22, Log-Lik: -1329.825, Max-Change: 0.02072
Iteration: 23, Log-Lik: -1329.815, Max-Change: 0.02112
Iteration: 24, Log-Lik: -1329.804, Max-Change: 0.02158
Iteration: 25, Log-Lik: -1329.734, Max-Change: 0.02587
Iteration: 26, Log-Lik: -1329.720, Max-Change: 0.02565
Iteration: 27, Log-Lik: -1329.707, Max-Change: 0.02592
Iteration: 28, Log-Lik: -1329.623, Max-Change: 0.03023
Iteration: 29, Log-Lik: -1329.607, Max-Change: 0.03013
Iteration: 30, Log-Lik: -1329.591, Max-Change: 0.03048
Iteration: 31, Log-Lik: -1329.494, Max-Change: 0.03517
Iteration: 32, Log-Lik: -1329.477, Max-Change: 0.03493
Iteration: 33, Log-Lik: -1329.459, Max-Change: 0.01506
Iteration: 34, Log-Lik: -1329.456, Max-Change: 0.04136
Iteration: 35, Log-Lik: -1329.437, Max-Change: 0.03616
Iteration: 36, Log-Lik: -1329.419, Max-Change: 0.03661
Iteration: 37, Log-Lik: -1329.309, Max-Change: 0.04138
Iteration: 38, Log-Lik: -1329.290, Max-Change: 0.04095
Iteration: 39, Log-Lik: -1329.271, Max-Change: 0.04120
Iteration: 40, Log-Lik: -1329.157, Max-Change: 0.04555
Iteration: 41, Log-Lik: -1329.137, Max-Change: 0.04525
Iteration: 42, Log-Lik: -1329.117, Max-Change: 0.04549
Iteration: 43, Log-Lik: -1329.002, Max-Change: 0.04943
Iteration: 44, Log-Lik: -1328.983, Max-Change: 0.04984
Iteration: 45, Log-Lik: -1328.964, Max-Change: 0.04835
Iteration: 46, Log-Lik: -1328.853, Max-Change: 0.05603
Iteration: 47, Log-Lik: -1328.834, Max-Change: 0.05338
Iteration: 48, Log-Lik: -1328.816, Max-Change: 0.05192
Iteration: 49, Log-Lik: -1328.712, Max-Change: 0.05430
Iteration: 50, Log-Lik: -1328.694, Max-Change: 0.05488
Iteration: 51, Log-Lik: -1328.677, Max-Change: 0.05551
Iteration: 52, Log-Lik: -1328.578, Max-Change: 0.05994
Iteration: 53, Log-Lik: -1328.562, Max-Change: 0.05838
Iteration: 54, Log-Lik: -1328.546, Max-Change: 0.05845
Iteration: 55, Log-Lik: -1328.456, Max-Change: 0.06069
Iteration: 56, Log-Lik: -1328.441, Max-Change: 0.06055
Iteration: 57, Log-Lik: -1328.427, Max-Change: 0.06088
Iteration: 58, Log-Lik: -1328.345, Max-Change: 0.06152
Iteration: 59, Log-Lik: -1328.332, Max-Change: 0.06212
Iteration: 60, Log-Lik: -1328.319, Max-Change: 0.06246
Iteration: 61, Log-Lik: -1328.245, Max-Change: 0.06481
Iteration: 62, Log-Lik: -1328.233, Max-Change: 0.06452
Iteration: 63, Log-Lik: -1328.222, Max-Change: 0.06456
Iteration: 64, Log-Lik: -1328.155, Max-Change: 0.06593
Iteration: 65, Log-Lik: -1328.145, Max-Change: 0.06591
Iteration: 66, Log-Lik: -1328.134, Max-Change: 0.06874
Iteration: 67, Log-Lik: -1328.073, Max-Change: 0.06981
Iteration: 68, Log-Lik: -1328.063, Max-Change: 0.06728
Iteration: 69, Log-Lik: -1328.054, Max-Change: 0.06611
Iteration: 70, Log-Lik: -1328.003, Max-Change: 0.07361
Iteration: 71, Log-Lik: -1327.994, Max-Change: 0.06790
Iteration: 72, Log-Lik: -1327.986, Max-Change: 0.06682
Iteration: 73, Log-Lik: -1327.940, Max-Change: 0.06446
Iteration: 74, Log-Lik: -1327.933, Max-Change: 0.07094
Iteration: 75, Log-Lik: -1327.926, Max-Change: 0.07007
Iteration: 76, Log-Lik: -1327.884, Max-Change: 0.07049
Iteration: 77, Log-Lik: -1327.877, Max-Change: 0.07142
Iteration: 78, Log-Lik: -1327.871, Max-Change: 0.08124
Iteration: 79, Log-Lik: -1327.831, Max-Change: 0.06560
Iteration: 80, Log-Lik: -1327.824, Max-Change: 0.06915
Iteration: 81, Log-Lik: -1327.818, Max-Change: 0.07046
Iteration: 82, Log-Lik: -1327.786, Max-Change: 0.07482
Iteration: 83, Log-Lik: -1327.780, Max-Change: 0.07337
Iteration: 84, Log-Lik: -1327.775, Max-Change: 0.07155
Iteration: 85, Log-Lik: -1327.746, Max-Change: 0.07769
Iteration: 86, Log-Lik: -1327.741, Max-Change: 0.07462
Iteration: 87, Log-Lik: -1327.736, Max-Change: 0.07441
Iteration: 88, Log-Lik: -1327.710, Max-Change: 0.07091
Iteration: 89, Log-Lik: -1327.706, Max-Change: 0.07165
Iteration: 90, Log-Lik: -1327.702, Max-Change: 0.07258
Iteration: 91, Log-Lik: -1327.678, Max-Change: 0.07570
Iteration: 92, Log-Lik: -1327.675, Max-Change: 0.07507
Iteration: 93, Log-Lik: -1327.671, Max-Change: 0.07409
Iteration: 94, Log-Lik: -1327.651, Max-Change: 0.07329
Iteration: 95, Log-Lik: -1327.647, Max-Change: 0.07495
Iteration: 96, Log-Lik: -1327.644, Max-Change: 0.07514
Iteration: 97, Log-Lik: -1327.625, Max-Change: 0.07639
Iteration: 98, Log-Lik: -1327.622, Max-Change: 0.07650
Iteration: 99, Log-Lik: -1327.619, Max-Change: 0.00312
Iteration: 100, Log-Lik: -1327.619, Max-Change: 0.07501
Iteration: 101, Log-Lik: -1327.616, Max-Change: 0.07486
Iteration: 102, Log-Lik: -1327.613, Max-Change: 0.07507
Iteration: 103, Log-Lik: -1327.597, Max-Change: 0.07353
Iteration: 104, Log-Lik: -1327.594, Max-Change: 0.07149
Iteration: 105, Log-Lik: -1327.592, Max-Change: 0.00414
Iteration: 106, Log-Lik: -1327.592, Max-Change: 0.00415
Iteration: 107, Log-Lik: -1327.591, Max-Change: 0.00021
Iteration: 108, Log-Lik: -1327.591, Max-Change: 0.07427
Iteration: 109, Log-Lik: -1327.589, Max-Change: 0.07333
Iteration: 110, Log-Lik: -1327.587, Max-Change: 0.07407
Iteration: 111, Log-Lik: -1327.584, Max-Change: 0.07417
Iteration: 112, Log-Lik: -1327.570, Max-Change: 0.07582
Iteration: 113, Log-Lik: -1327.568, Max-Change: 0.07539
Iteration: 114, Log-Lik: -1327.566, Max-Change: 0.05673
Iteration: 115, Log-Lik: -1327.561, Max-Change: 0.00505
Iteration: 116, Log-Lik: -1327.560, Max-Change: 0.08820
Iteration: 117, Log-Lik: -1327.558, Max-Change: 0.07629
Iteration: 118, Log-Lik: -1327.546, Max-Change: 0.00593
Iteration: 119, Log-Lik: -1327.545, Max-Change: 0.00122
Iteration: 120, Log-Lik: -1327.545, Max-Change: 0.07510
Iteration: 121, Log-Lik: -1327.544, Max-Change: 0.00080
Iteration: 122, Log-Lik: -1327.544, Max-Change: 0.07048
Iteration: 123, Log-Lik: -1327.542, Max-Change: 0.06928
Iteration: 124, Log-Lik: -1327.532, Max-Change: 0.06730
Iteration: 125, Log-Lik: -1327.530, Max-Change: 0.00117
Iteration: 126, Log-Lik: -1327.530, Max-Change: 0.00300
Iteration: 127, Log-Lik: -1327.530, Max-Change: 0.00188
Iteration: 128, Log-Lik: -1327.530, Max-Change: 0.00069
Iteration: 129, Log-Lik: -1327.530, Max-Change: 0.00076
Iteration: 130, Log-Lik: -1327.530, Max-Change: 0.00018
Iteration: 131, Log-Lik: -1327.530, Max-Change: 0.00055
Iteration: 132, Log-Lik: -1327.530, Max-Change: 0.00057
Iteration: 133, Log-Lik: -1327.530, Max-Change: 0.00023
Iteration: 134, Log-Lik: -1327.530, Max-Change: 0.00048
Iteration: 135, Log-Lik: -1327.530, Max-Change: 0.00014
Iteration: 136, Log-Lik: -1327.530, Max-Change: 0.00011
Iteration: 137, Log-Lik: -1327.530, Max-Change: 0.00251
Iteration: 138, Log-Lik: -1327.530, Max-Change: 0.00075
Iteration: 139, Log-Lik: -1327.530, Max-Change: 0.00094
Iteration: 140, Log-Lik: -1327.530, Max-Change: 0.00015
Iteration: 141, Log-Lik: -1327.530, Max-Change: 0.00045
Iteration: 142, Log-Lik: -1327.530, Max-Change: 0.00058
Iteration: 143, Log-Lik: -1327.530, Max-Change: 0.00016
Iteration: 144, Log-Lik: -1327.530, Max-Change: 0.00043
Iteration: 145, Log-Lik: -1327.530, Max-Change: 0.00093
Iteration: 146, Log-Lik: -1327.530, Max-Change: 0.00054
Iteration: 147, Log-Lik: -1327.530, Max-Change: 0.00027
Iteration: 148, Log-Lik: -1327.530, Max-Change: 0.00022
Iteration: 149, Log-Lik: -1327.530, Max-Change: 0.00038
Iteration: 150, Log-Lik: -1327.530, Max-Change: 0.00013
Iteration: 151, Log-Lik: -1327.530, Max-Change: 0.00011
Iteration: 152, Log-Lik: -1327.530, Max-Change: 0.00041
Iteration: 153, Log-Lik: -1327.530, Max-Change: 0.00038
Iteration: 154, Log-Lik: -1327.530, Max-Change: 0.00021
Iteration: 155, Log-Lik: -1327.530, Max-Change: 0.00041
Iteration: 156, Log-Lik: -1327.530, Max-Change: 0.00013
Iteration: 157, Log-Lik: -1327.530, Max-Change: 0.00011
Iteration: 158, Log-Lik: -1327.530, Max-Change: 0.00188
Iteration: 159, Log-Lik: -1327.530, Max-Change: 0.00071
Iteration: 160, Log-Lik: -1327.530, Max-Change: 0.00278
Iteration: 161, Log-Lik: -1327.530, Max-Change: 0.00053
Iteration: 162, Log-Lik: -1327.530, Max-Change: 0.00011
Iteration: 163, Log-Lik: -1327.530, Max-Change: 0.00009
##
## Calculating information matrix...
m_depression
##
## Call:
## mirt(data = df_dep, model = 1, itemtype = "graded", SE = TRUE)
##
## Full-information item factor analysis with 1 factor(s).
## Converged within 1e-04 tolerance after 163 EM iterations.
## mirt version: 1.37.1
## M-step optimizer: BFGS
## EM acceleration: Ramsay
## Number of rectangular quadrature: 61
## Latent density type: Gaussian
##
## Information matrix estimated with method: Oakes
## Second-order test: model is a possible local maximum
## Condition number of information matrix = 72063
##
## Log-likelihood = -1328
## Estimated parameters: 10
## AIC = 2675
## BIC = 2717; SABIC = 2685
## G2 (14) = 27.1, p = 0.0188
## RMSEA = 0.044, CFI = NaN, TLI = NaN
# Factor loadings
summary(m_depression)
## F1 h2
## Depression_3 0.995 0.989
## Depression_4 0.817 0.667
##
## SS loadings: 1.66
## Proportion Var: 0.828
##
## Factor correlations:
##
## F1
## F1 1
# Alpha
a2 <- coef(m_depression)
a2 <- data.frame(
Item = c("Depression 1", "Depression 2"),
Discrimination = c(paste0(format_value(a2$Depression_3[1, 1]), ", ", format_ci(a2$Depression_3[2, 1], a2$Depression_3[3, 1])), paste0(format_value(a2$Depression_4[1, 1]), ", ", format_ci(a2$Depression_4[2, 1], a2$Depression_4[3, 1]))
))
knitr::kable(a2)
| Item | Discrimination |
|---|---|
| Depression 1 | 16.46, 95% CI [-10.08, 43.00] |
| Depression 2 | 2.41, 95% CI [2.03, 2.79] |
out <- crc(m_depression, df_dep)
p2b <- out$p + labs(x = expression(Depression~(theta)))
p2b
p3b <- normalize_scores(out$rez, out$minmax) + labs(x = expression(Depression~(theta)))
p3b
fig1 <- wrap_elements(p1) + wrap_elements(p2a / p2b + plot_annotation(title = "B. Information Curves per Response Type", theme = list(plot.title = element_text(face = "bold"))) + plot_layout(guides = "collect")) + plot_layout(widths = c(1, 1.5))
ggsave("figures/figure1.png", fig1, width = fig.width*1.43, height=fig.height*1.43)
fig2 <- p3a / p3b + plot_annotation(title = "Normalized Scoring", theme = list(plot.title = element_text(face = "bold"))) + plot_layout(guides = "collect")
ggsave("figures/figure2.png", fig2, width = fig.width*1.43, height=fig.height*1.43)
data <- sapply(cbind(df_anx, df_dep), function(x) {
case_when(
x == 0 ~ "Not at all",
x == 1 ~ "Once or twice",
x == 2 ~ "Several days",
x == 3 ~ "More than half the days",
TRUE ~ "Nearly every day"
)
}) |>
as.data.frame()
source("score_PHQ4.R")
data <- score_PHQ4(A1=data$Anxiety_1, A2=data$Anxiety_2, D1=data$Depression_3, D2=data$Depression_4, method="basic") |>
datawizard::data_addsuffix("_Basic") |>
cbind(
score_PHQ4(A1=data$Anxiety_1, A2=data$Anxiety_2, D1=data$Depression_3, D2=data$Depression_4, method="normalized") |>
datawizard::data_addsuffix("_Normalized")) |>
select(matches("Depression|Anxiety")) |>
pivot_longer(everything()) |>
separate(name, into = c("Dimension", "Scoring"))
data |>
ggplot(aes(x=value)) +
geom_area(data=estimate_density(data, at = c("Dimension", "Scoring"), method = "KernSmooth") |>
group_by(Scoring) |>
normalize(select = "y"),
aes(x = x, y = y, fill=Dimension), alpha = 0.05, position="identity") +
geom_hline(yintercept = 0.5, linetype = "dotted") +
stat_ecdf(aes(color=Dimension), geom = "step", size=1) +
facet_wrap(~Scoring, scales = "free") +
scale_y_continuous(labels=scales::percent, expand = c(0, 0)) +
scale_x_continuous(expand = c(0, 0)) +
scale_fill_manual(values = c("Depression" = "#3F51B5", "Anxiety" = "#E91E63"), guide = "none") +
scale_color_manual(values = c("Depression" = "#3F51B5", "Anxiety" = "#E91E63")) +
theme_modern(axis.title.space = 5) +
labs(title = "Distribution and Cumulative Density of the Sample", x = "Score", y = "Cumulative Proportion")
df$Depression <- mirt::fscores(m_depression)[, 1]
df$Anxiety <- mirt::fscores(m_anxiety)[, 1]
correlation::correlation(df) |>
summary() |>
export_table(caption = "Correlation Matrix")
## Correlation Matrix
##
## Parameter | Anxiety | Depression | Depression_4 | Depression_3 | Anxiety_2
## -----------------------------------------------------------------------------
## Anxiety_1 | 0.89 | 0.72 | 0.56 | 0.72 | 0.83
## Anxiety_2 | 0.98 | 0.72 | 0.57 | 0.73 |
## Depression_3 | 0.74 | 0.98 | 0.73 | |
## Depression_4 | 0.59 | 0.78 | | |
## Depression | 0.75 | | | |